A DIRICHLET-TO-NEUMANN APPROACH FOR THE EXACT COMPUTATION 
OF GUIDED MODES IN PHOTONIC CRYSTAL WAVEGUIDES 



SONIA FLISSt 

Abstract. This works deals with one dimensional infinite perturbation - namely line defects - in periodic media. 
In optics, such defects are created to construct an (open) waveguide that concentrates light. The existence and the 
computation of the eigenmodes is a crucial issue. This is related to a self-adjoint eigenvalue problem associated 
to a PDE in an unbounded domain (in the directions orthogonal to the line defect), which makes both the analysis 
and the computations more complex. Using a Dirichlet-to-Neumann (DtN) approach, we show that this problem is 
equivalent to one set on a small neighborhood of the defect. On contrary to existing methods, this one is exact but 
there is a price to be paid : the reduction of the problem leads to a nonlinear eigenvalue problem of a fixed point 
nature. 
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1. Introduction. Periodic media play a major role in applications, in particular in op- 
tics for micro and nano-technology ll20l ISTl l27l [38l . From the point of view of applications, 
one of the main interesting features is that it can exist intervals of frequencies for which the 
propagative waves cannot exist in the media. This phenomenon is due to the fact that a wave 
on the media is multiply scattered by the periodic structure, which can lead, depending on the 
characteristics of the media and the frequency, to possibly destructive interferences. It seems 
then that an appropriate choice of the structure and the dielectric materials of the photonic 
crystal can create particular band gap and then, from a practical point of view, banish some 
monochromatic electromagnetic waves. Thus, the periodic media could be used to several 
potential applications such as in the realization of filters, antennas and more generally, com- 
ponents used in telecommunications. 

Mathematically, this property is linked to the gap structure of the spectrum of the under- 
lying differential operator appearing in the model. For a complete, mathematically oriented 
presentation, we refer the reader to [27 , 28 1. Even if the necessary conditions for the existence 
of band gaps are not known -except for the one dimensional case in [4J where authors state 
that the absence of gaps implies that the coefficients of the media are constant-, sufficient 
conditions exist. Figotin and Kuchment have given examples of high contrast medium for 
which at least one band gap exists and can be characterized lfT3l[T4l . Using asymptotic argu- 
ments, Nazarov and co-workers have established that a small perturbation of an homogeneous 
waveguide can open a gap in the continuous spectrum of the operator l.31i i7J. Moreover, other 
band gap structures have been characterized through numerical approaches in jS). 

Besides, Sommerfeld and Bethe in 1933 made the conjecture that the number of gaps has 
to be finite in 2d and 3D. This conjecture has been proven by Skriganov for Schrodinger op- 
erator in 2D BOl [36J and in 3D BTI l42l - the last article goes further: if the contrast is too 
small, gaps cannot exist and on suitable assumptions on the lattice, the conjecture holds for 
higher dimensions. Parnovski |32| has generalized this last result provided that the potential 
is sufficiently smooth and in Il33l the conjecture was proven in full generality. Another proof 
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can be found in the book of Yu Karpeshina ll24l . For Maxwell equations, it is much more 
complicated. Nevertheless, there is now evidence that gaps exists. Experimentally this was 
first observed by Yablonovitch et al |47l. There are also many numerical results for different 
structures (see |43, 44] for references). Recently, Vorobets in |46| has proven the conjec- 
ture for 2D periodic maxwell operator (which corresponds to a 2D periodic photonic crystal 
which is homogeneous in the third direction) with separable dielectric function. Moreover, if 
the dielectric function is close enough to a constant, there is no gaps at all. 

These media can present perturbations or defects which are introduced in the media to change 
their properties. Thus, in optic, in order to produce lasers, fibers or waveguides in general, 
it is necessary to have authorized frequencies inside the forbidden intervals of frequencies. 
This property can be obtained, for example, introducing localized defect or line defect and 
corresponds, from a mathematical point of view, to isolated eigenvalues of finite multipUcity 
inside the gaps of the underlying differential operator appearing in the model. 

Figotin and Klein have proven rigorously that introducing a defect in a periodic structure, 
namely a perturbation of compact support, can create defect modes, which are eigenvectors 
associated to eigenvalues inside a gap lfT0l[Tn[T2l . More precisely, these defect modes are 
stationary waves exponentially decreasing far from the perturbation, then they seems to "live" 
around the defect which explains their name. Using asymptotic arguments, Nazarov has stud- 
ied in 1 30 1 sufficient conditions for existence of eigenvalues below the essential spectrum for 
elastic waveguide. Here again, it seems difficult to obtain necessary conditions for the exis- 
tence of defect modes. 

There are very few works in the same spirit in the mathematical literature for the case of 
line defect, that means a one dimensional perturbation of the periodic medium. A first natural 
question concerns the spectrum of the perturbed operator compared to the spectrum of the 
perfectly periodic operator. More precisely, if a gap is in the spectrum of the perfectly peri- 
odic operator, can a part of the spectrum of the perturbed one arise in this gap? Moreover, if a 
part of a spectrum arises, is it of absolutely continuous or singular type? Does it corresponds 
to guided modes, meaning that the waves are confined to the guide (evanescent in the periodic 
media) and they are propagating along the line defect? Kuchment et Ong 1291 have answered 
partially to these questions for the case of homogeneous line defect. Indeed, for any gap (a, b) 
of the perfectly periodic operator, they found a sufficient condition on a, b, the characteristics 
of the line defect (its size and its coefficients) such that spectrum of the perturbed operator 
arises in the gap. Moreover, they have shown that the associated eigenvector is confined to 
the guide. However, they have not concluded on the existence or not of bound states - this 
question is linked to the nature of the spectrum. Vu Hoang and Maria Radosz in |19| have 
shown for waveguide in 2D periodic structures the absence of bound states. It seems then that 
if the spectrum arises in the gap, it corresponds to guided modes. Let us remark that for the 
particular case of a line periodic perturbation in a homogeneous media, Bonnet-BenDhia and 
Starling |3 1 have given necessary and sufficient conditions for the existence of guided modes. 

Finally, few works (| 1 1 and the present one) characterized precisely (respectively using Green's 
function or DtN operator) the guided waves if they exist, for general line defect. These par- 
ticular modes can help for the determination of the spectrum of the perturbed operator. The 
advantage is that the study can be reduced to a band (bounded in the direction of the period- 
icity and infinite in the other direction). 
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From a numerical point of view, there exist only few methods. The most known is the Su- 
percell method. It consists in making computations in a bounded domain of large size with 
periodic boundary conditions, the resulting solution converging to the true solution when the 
size tends to infinity. The convergence for the computation of defect modes has been shown 
for 2D problems and compact perturbations in ll45l and generalized to 3D problems and to 
exponentially decreasing perturbations in [6'|. In this case, as the localized modes are ex- 
ponentially decreasing, this convergence is exponentially fast with respect to the size of the 
truncated domain. In practice, this approach replaces the eigenvalue problem set in an un- 
bounded domain to an approximated one set in a bounded domain. See |39 | for numerical 
results. The main drawback of this strategy relies on the increase of the computational cost, 
especially when a mode is not well confined. We can mention also the fictitious source su- 
perposition method fSl and the reflective scattering matrix method. 

By adapting to eigenvalue problems the construction of Dirichlet-to-Neumann operators orig- 
inally developed for scattering problems ll22l[T6l[T5]| . we want to offer a rigorously justified 
alternative to existing methods. Compared to the supercell method, the DtN method allows 
us to reduce the numerical computation to a small neighborhood of the defect independently 
from the confinement of the computed guided modes. Moreover, as the method is exact, we 
improve the accuracy for non well-confined guided modes. Obviously, there is a price to be 
paid : the reduction of the problem leads to a non linear eigenvalue problem, of a fixed point 
nature. However, this difficulty has been akeady overcome for homogeneous open waveg- 
uides for which the DtN approach is well known ll23l [34l |35]| . 

2. Model problem. In order to describe the medium of study, let us first consider a 
two-dimensional periodic medium - the photonic crystal - characterized by a coefficient pp 
(typically the square of the refraction index of the medium) which is 

• a. L°° function : 3p-,p+, < p_ < pp{x,y) < p+; 

• periodic in the two directions 

3Lx,Ly > 0, Vn, m e Z, V(x, y) e pp{x + nL^^y + mLy) = Pp{x,y) 

where the periods and Ly are not necessarily equal. We introduce a one dimensional infi- 
nite perturbation - called the line defect - in the y— direction in 

flo = ] — a, a[ X M 

and characterized by a coefficient po which is 

• a L°° function satisfying < ^ Po 1^ P+'^ 

• periodic in the y— direction 

3Ll > 0, Vm e Z, yix, y) e po{x, y + mil) - pa{x, y) 
where Ly and are commensurate. Without lack of generality, we suppose Ly = 

Ly. 

The propagation medium is then characterized by the function p defined by 



p{x,y) 



Pp{x,y) in M- C f^o 
Po{x,y) in ■ 



Remark 2.1 (Some extensions). 
L We could consider a more general medium of propagation than fl = M^. It could consist 
for example o/M^ minus a periodic set of holes. The only assumption is that the periodicity 
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Fig. 2.1. Domain of propagation : typically p = 1 in the white region, p = 2 in the dark grey regions and 
p = 3 in the light grey regions. 

property of^ has to the be the same as p. 

2. Even if the spectral properties are different, the method can apply to coefficient defined by 



where and have the same properties than pp with the same period in the y— direction 
(Ly) and not necessarily the same in the x-direction ( and ). 

Remark 2.2 (Some open questions). The case where the line defect has not the same 
periodicity properties than the photonic crystal and the case where the line defect is not in- 
troduced along one of the direction of periodicity of the photonic crystal are out of the scope 
of the present paper. 

The propagation model is a simple 2D space-(x = (x, y)) time harmonic scalar wave equa- 
tion (corresponding for example to transverse electric (TE) mode electromagnetic wave in 
two dimensions) 



where u G M+ is the frequency. In the rest of the article, uj'^ will be considered as the spectral 
parameter. 

Remark 2.3. The results developed in this article can be easily extended to more gen- 
eral elliptic operator u i— > V • (/iVu) where /i is line perturbation of a periodic function. 

Using the Floquet-Bloch theory, we could show that the spectrum of the operator — A/p 
could be deduced from the study of the guided modes. A guided mode of this problem is by 
definition a solution w 7^ to which can be written in the form 



Pp{x,y) in ]— 00,— a[x]R 
Po{x,y) in no 
Pp(x,y) in ]a, +oo[xIR. 



Aw — pJ^ w = Q, in il, 



w{x,y) = v{x,y) e'' 



(2.2) 
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Ly Ly 



where, in full generality 

• /3 - called the quasi period - is in ] — n/ Ly^n/ Ly[; 

• w is periodic in the direction with period Ly and 

v\^eH\B) where B = Mx; 2' 2' 

We will identify in the following any periodic function in il to its restriction to the period B. 

Remark 2.4. In full generality, the coefficient /? is in M but it is enough to consider it 
only in ] — / Ly, tt / Ly[. Indeed, it is easy to see that if there exists a mode for /3 G M, 

'w{x, y, t) = v{x, y) e^^^ , with v Ly-periodic and lu G M^, 

it corresponds to a mode for /3 ± /Ly of the form 

w{x, y) = v{x, y) e^^^^'^/^y^'^ with v{x, y) = v{x, y)e^^'^/^^y which is Ly-periodic. 

Replacing the expression \2.2\ in (|^ we see easily that finding the guided modes corresponds 

r- r+ y. 



B. 



O O OiC^iO O Ob, 



Fig. 2.2. The hand B - one period of the domain in the y-direction. 

to finding couples (o;^, /3) such that there exists v E H^{B), v ^ solution of 

1 



-Av ■ 



P 



2iP ^2 

vv H V = uj V, in B 

P P 



dyV 



(2.3) 



where (see figure 2.2 1 E = M x {^} and S = M x { - ^^}. 



We can rewrite the previous problem using u{x, y) = v{x, y)e^^y and then finding the guided 
modes is equivalent to finding couples (w^,/?) such that there exists u E H^(B), u 7^ 
solution of 



— -Au = oj'^u, in B 
P 



dyU I5. 



= e"^ a, u 



(2.4) 



To characterize the guided modes, there exist two different formulations : the w— formulation 
in which cj is fixed and /3 is looked for and the /3— formulation in which f3 is fixed and u is 
looked for. The first one leads to a quadratic eigenvalue problem in the band B (it is obvious 
with the formulation ( |2.3| ) of the problem) and the second one to an eigenvalue problem in 
the band B. 



The classical approach for solving this problem set in an unbounded domain is the Super- 
cell Method which is based on the exponentially decreasing property of the mode in the x— 



6 



S. Fliss 



direction (see Theorem 3.5 where we remind this result) and consists in truncating the band 
B far enough. One (artificial) parameter of the method is then the position of the truncature 
position. The main advantage is that this approach leads to an eigenvalue problem (for the 
/3-formulation) or a quadratic one (for the cj— formulation) set on a bounded domain (the trun- 
cated one). The main drawback of this strategy relies on the increase of the computational 
cost when a mode is not well confined. 



Here we propose a novel method based on a DtN approach, which is originally developed for 
scattering problems. This method offers a rigorously justified alternative to existing meth- 
ods. This approach allows us to reduce the numerical computation to a small neighborhood of 
the defect independently from the confinement of the computed guided mode. The main ad- 
vantage is that the method is exact so we improve the accuracy for non well-confined guided 
modes and could deduce the behavior of the dispersive curves - (3 where {(3,uj{/3)) 

is solution of ( |2.3| l - near the edges of the gaps. The main drawback is that the reduction of 
the problem leads to a nonlinear eigenvalue problem, of a fixed point nature. To simplify the 
presentation, we choose here the formulation but the method extends to the other formu- 
lation - which could be more adapted for dispersive media for example, p — p{i^)). 

The article is organized as follows. We will remind in Section 3 the spectral properties of 
the locally perturbed periodic operator involved in the study. Section 4 deals with the DtN 
approach and present the non linear eigenvalue problem which has to be solved to compute 
the guided modes. It is the most important and original part of the article. Section 5 is devoted 
to numerical results and Section 6 to give some conclusions and perspectives. 

3. Spectral theory results for the problem. We focus on the guided modes defined in 
\2A\ and then reduce the problem to the following one : 

Forany ^ e] - 

l^y l^y 

findw2£M+, s.t. 3u e u 7^ 0, A{l3)u = uj'^u (£^) 

where 

A{/3)^--A 
P 

D{A{f3))^{ueH^{A,B), w|s=e''^u|g and dyu\^ = e'^dyu\^}. 
with H^iA, B) = {ue H^iB), Au e L^iB)}. 



The problem is then reduced to the determination of eigenvalues of the locally perturbed 
periodic operator A{f3), for any /3 e] — tt/ Ly,TT / Ly[. 

We introduce for the following the operator with perfectly periodic coefficients 

Ap{(3) = --A, D{Apm=D{Am (3.1) 
Pp 

Using fSS^ TOl and the Weyl's theorem we describe in the following proposition the essential 
spectrum of the operator A{/3). 
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Proposition 3.1 (Essential spectrum of A{f3)). The operator A{j3) is selfadjoint in 
H = L'^{B, pdxdy), positive and its essential spectrum, denoted (Jess{P)t satisfies 

aess{P)^<j{Ap{(i))^R\ y ]a„(/3),6„(/3)[ 

neIl.Af(/3)l 

where Ap{j3) is the operator with periodic coefficient defined by p.l| l, and N{j3) (0 < 
< +ooj is the number of gaps, < an{(3) < bn{l3). 

Proof. The resolvant operator of A{i3) is a compact perturbation of the one of Ap{f3), so 
the essential spectrum of A{/3) coincide with the essential spectrum of Ap{f3) 



(3.2) 



Moreover, one of the main result of the Floquet-Bloch theory is (see 11261 for more details) 
that the spectrum of Ap{p), a {Ap{P)) is reduced to its essential spectrum, ciess (^p(/3)) and 
is given by 



<j{Apm^<j,,,{Apm^ y <T{Ap{(3,k)) 



(3.3) 



where 



Ap{p,k) = A, 

Pp 



D{Ap{p,k)) = {ueH\A,C), 



" \x. 

dxU 
u ' 



— ptkL^ I 



\x=L^/2 - "2^" \x=-Lj2 

_ = " " \y=-Ly/2 ' 

withC =] - Lx/2,Lx/2[x] - Ly/2,Ly/2[iindH\A,C) ^ {u e H^{C), Au e L^{C}. 

Moreover, for any fc in ] — tt/Lj,, tt/L^], Ap{(5, k) is a self-adjoint positive operator with 
compact resolvent so its spectrum is purely discrete 

< wi(/3, fc) < a;2(/9, fc) < . . . < w„(/3, fc) < . . . with lim a;„(/3, fc) = +cx) 

and we can find (e„(/3, k)) a corresponding family of eigenvectors which is an hilbertian ba- 
sis of L2(c). 

Using the min-max principle to Ap{f3, k), we could show that the so called dispersive curve 

k) I— 7> uJn{(3, k) is continuous. (3.4) 
Finally, we deduce from p.2| ), \33) and p.4| i that 

a,,,(/3)=a(Ap(/3))= y c^„([-^, /?) 



□ 

Remark 3.2. In Section 4.4 we will give another characterization of the essential spec- 
trum with a by-product of the method. 
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As described in the Introduction, only sufficient conditions of existence are known and finite- 
ness of the number of gaps is conjectured in the general case and proven for particular cases. 
Let us suppose in the following that at least one gap exists iN{/3) < 1). Using the theory of 
selfadjoint operators ( ll37]| ). we deduce that 

Proposition 3.3 (Discrete spectrum of A{/3)). The spectrum ofA{(3) inside the gaps 
consists only of isolated eigenvalues of finite multiplicity, which can accumulate only at the 
edges of the gap. 

One should ask if there is a way to ensure the rise ot at least one eigenvalue in gaps of 
Figotin and Klein have given in |10| sufficient condition on the defect (supposed ho- 
mogeneous) and the gaps to introduce eigenvalues in any gaps of A(/3). Using asymptotic 
arguments, Nazarov and co-workers 1301 have studied sufficient conditions for existence of 
eigenvalue for elastic waveguides. 

We are interested now, in characterizing and then computing the eigenvalues (Am(/?))m, 
< Ai (/?) < A2 (/3) < . . . < Am(^,) (/3) , < M{P) < +(x) 

if they exist (we suppose then Af(/3) > 1), which are in the gaps of the essential spectrum (see 
I -optical waveguides- and ||29J -photonic crystal waveguides- for existence of eigenvalues 
inside gaps) : 

Proposition 3.4 (Properties of each A™ (/?)). The dispersion curves (3 i-> \m{(i) are 
2t: / Ly-periodic, even and continuous. 

Proof. 

• /? A„i(/?) is 27r/Lj^ -periodic : by definition 

A{f3 + 271 /Ly)=A{P) and D (A{P + 2Tr/Ly)) = D {A{/3)) . 

• /? I— > A„j(/?) is even because 

Ail3)=A{-P) and D {Am = D {A{-p)) . 

and A{/3) is a self-adjoint positive operator. 

• The continuity with respect to /3 of the eigenvalues sorted in ascended order is due 
to the analyticity of the operator with respect to (3 ||25). 

□ 

We deduce in particular that it is sufficient to study the dispersive curves for f3 E [0, n/Ly]. 

We study now the properties of the eigenvectors associated to the eigenvalues. By defini- 
tion, they are in L^{B) but we can be more precise : they decay exponentially fast far from 
the defect, with a rate depending on the distance from the eigenvalue to the edges of the gap. 
It is a result based on ||9l and shown in lITOl . 

Theorem 3.5 (Exponential decay of the eigenvectors). Let /3 e] — n / Ly^n / Ly\. For 
any eigenvalue A„(/3) of A(/3), the associated eigenvectors ipniP', •) satisfies 

Vn e Il,M(/3)], 3Ci,C2 > 0, 3a > 0,V(a;,zj) G B, \x\ < a, 

^ r/' e(-^--^"'(^-(^)---(^(^)))l-l) (3.5) 

dist{Xn{l3),aess{A[fi))) 
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This propriety is exactly the one which encourages to use the Supercell method, hideed, 
this method consists in approaching the eigenvalues and the associated eigenvectors by the 
ones of a truncated problem in the x— direction with periodic conditions. If the eigenvec- 
tors are exponentially decreasing in the x— direction, they will satisfy almost periodic condi- 
tions far enough from the perturbation and then be a/mos?-eigenvectors for the almost-same 
eigenvalues of an operator defined from the truncated domain. The main advantage of this 
approach it that it leads to eigenvalue problem for a /3-formulation and a quadratic one for 
a oj— formulation, both of them set on a bounded domain. However, we could note several 
drawbacks. 

1. The essential spectrum of A{f3) has to be computed initially. 

2. For a fixed /3, it seems important to have an estimation of the distance between the 
not yet computed eigenvalue and the essential spectrum of A{(3) to choose a relevant 
truncated domain. 

3. The size of the truncated domain depends on /3 and on the distance between the 
eigenvalue and the essential spectrum. If the eigenvalue approaches more and more 
the essential spectrum of A{f3) when /3 varies, the corresponding eigenvector be- 
comes less and less confined (less and less exponentially decreasing) and then the 
truncated domain has to be bigger and bigger This can increase dramatically the 
computational cost. 

4. Because, this method is based on the exponential decay of the eigenvectors, it cannot 
describe the behavior of the dispersive curves (the eigenvalues as functions of f3, 
P I— >■ A„(/3)) near the edges of the essential spectrum. 

We propose now a novel method based on a DtN approach which overcome these disad- 
vantages because it is an exact method (in the sense which is precised in the next section). 
However, its main drawback is that it leads to a nonlinear eigenvalue problem of a fixed point 
nature. 

4. The non linear eigenvalue problem. 

4.1. The DtN approach. Let /3 g] — tt / Ly^ir / Ly] and let us suppose from this point 
that ^ a^ssiP)- We want here to write an equivalent problem to the problem {Ep i which 
is set on a bounded domain. 



For the sequel, it is essential to introduce functional S£aces appearing naturally in the study. 
Let us remind that B± = S n f]=^, E± = S n 1]=^, E± = E n and r± = {±a}x] - 
iy/2, Ly/2[ (see figure 2.2 1. 



We start from smooth quasi-periodic functions m B^: 

i|^±, C°°(r!±), u{x,y + L) 

Let Hl{B^) be the closure ofCf{B^) in H^{B^). 

Hl{B^) = [u e H\B^), u\^^ = e^P^y 

where in the last equality we have identified the spaces iJ^/^(E*) and i7^/^(E*). As 
H^{B^) is a closed subspace of H^{B^), we equip it with the norm of H^{B^). Let 
i7^(A,B±)betheclosureof i7i(A,B±) = {m e H^{B^), Au e L2(S±)}. 



du 
dy 



du 
dy 



}■ 
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where in the last equality we have identified the spaces HqI^'^{Y,^)' and Hqq'^{Y,^)' . 
The space i7y^(r±)is defined by 

where 7^ € i/i/2(r±)) is the trace map on r± : Vu G H^{B^), IqU = u|p±. 

H^J^{11q) is then a dense subspace of _ff^/^(So) and the injection from i/^^^(rj) onto 
JJ^/^(r^) is continuous. 

We define H^^''^{T^) as the dual ofH]/'^{T^). 

Finally, the trace application 7^ e (A, i/^^a) (rj )') defined by : 

V,ieHi(A,S±), 7^=^ ± 
is a continuous application from iJ^(A, _B*) onto i/^ ^^'^(rj) and we can show that 

V^'(rJ)=7i(ffMA,i3±)). 
Let us now introduce the two half-band problems: for any /? and a and for any given ip in 



Findu± e H}i{A,B^), 



-Am± - p„q;2u=^ = in B± 
m|p± = (p. 



Theorem 4. 1 (Well-posedness of the problem ( |7^ 

Ifa^ ^ aess{l3), the problem {V^ \ is well-posed in Hjj except for a countable set of frequen- 
cies which depends on (3. 

If the periodicity cell is symmetric with respect to the axis x = and if ^ dess (/?), 
the problem (V^ I is always well-posed in H^. 

Proof. It is enough to show the result for (V^). We introduce the following operators 

Ad,+ W) = -^A, DiAD,+ W)) = e H'p{B+), u+\^+ = 0} 
Ad,-{P)^-^A, DiAD,-m^{u- eH}j{n\B+), M-|p+-o} 
AdW) = "^A, DiAom = {ue Hlin), u\^+ = 0} 

It is easy to see that the spectrum of Ajj is the union of the spectrums of and _. 

Moreover the resolvant of is a compact perturbation of the resolvant of A so and 
A has the same essential spectrum. These two properties implies that + has its essen- 
tial spectrum included in the spectrum of A but this operator may have a countable set of 
eigenvalues. In conclusion, ^ (Tess(/?), the problem (V^) is well posed if is not an 
eigenvalue of 
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If the periodicity cell is symmetric with respect to the axis a: = 0, the operator has 
no eigenvalues. Indeed, if an eigenvector exists, we could construct by antisymmetry with 
respect to r+ an eigenvector of A, which is not possible. □ 



Remark 4.2 (Robin-to-Robin operators instead of DtN operators). The countable set of 
frequencies for which the half-band problems are not well-posed is introduced by our method 
of construction of DtN operators. If we use Robin-to-Robin operators instead of Dirichlet-to- 
Neumann operators, the problems are always well posed when the frequency c? ^ cress(/3)- 
We choose to use DtN operators in the present paper to simplify the presentation of the method 
(see 47511781/ for more details on Robin-to-Robin operators) 

Suppose from this point that the problems {V^) and [V^) are well posed. We denote by 
it+(/3, a; ip) and (/3, a; (f) the respective unique solutions. 

The DtN operators A±(^, a) e C{hI''^{T^),H^^''^{T^)) are then defined by 



orelseV^,VG77y'(rJ) 



< A±(^,a)(/3,V >r±= / Vu*(/3,a;v3) • Vu*(/3,a;?/') 



-aM ppU±(/3,a;(^)u±(/3,a;V)- (4.1) 
Jb± 

where < •, • >p± denotes the duality product between ^^^(Fj) and H^J'^iV^). We have 
then by definition the continuity properties of the DtN operators 

Proposition 4.3. For any j3, the DtN operators K'^{P, a) are continuous from H^^'^{T'^) 

onto Hp ^^^(rj ) and are norm continuous with respect to a. 
Proof. We introduce the operator 

AdAP) = -^A, D{AdAP)) = {^^ e 7Ji(73±), u\i_ = o} 

We have seen in Theorem 



4.1 



that if ^ aessiP) and ^ (^d,±{P) where cFd.±{P) is the 
discrete spectrum of AuAfi) then the problem (P*) is well posed. Moreover, it exists C 
independent from a such that 

C 

V(pe i7y'(rj), ||K±(/?,a;^)||Hi(B±) < ^-p-^^^||<p||^i/2(j,±) (4.2) 



and so 

C 

y^^H}J\T% ||A±(/3,a)^||^-v.(p±) < dist(«2,,^(^))ll^ll<^r^)- 
where a± {/3) is the spectrum of Ad.± (/?). 
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Let now ai and a2 are such that each corresponding half-band problem (V^) is well posed. 

For any /3 and any data ip in hI^^{T^), — (/3, ai;ip) — u^(/3, 012; (p) is solution of 



- Ppofv^ — Pp{a2 ~ af)u^ {(3 , a2', f) in i?^ 
0. 



That implies that it exists C independent from ai and a2 such that 

C 

"""^^("")-distK,a±(/3))'"'" 
then using ( |4.2| i, we deduce that V(p G //^^(rj), 



a^|||M=^(^,a2;(/?)||m(B±) 



|A±(/3,ai)^-A±(/3,a2)^||^ 



< 



C 

.2 ^. ^« ^^^;..^^,2 ^. 



dist(Q;J, CT± (/?)) dist(a^, cr±(/3)) 



□ 

Using the continuity of any eigenvector of ( £^ 1 and its normal derivative across the sections 
r+ and r~, the next theorem is therefore straightforward. 



Theorem 4.4 (Problem with DtN conditions). The problem iSjS 1 is equivalent to the 
problem posed on Bq = B n CIq 

Find ui^ ^ (Tf,gg{(3), s.t. Buq £ H^{Bq), uq ^ 



-Alio = ^"^uq, in Bq 
P 



(^o) 



uq satisfying the boundary conditions 

+9^Mo + A+(/3,w)uo = 0, onV+ 
-(9a;Wo + A"(/3,a;)uo = 0, onT~, 

where Eo = S H f^o and Sq = S n JIq- These problems are equivalent in the sense that if 
(w, u) is solution of 1 then (oj, u\bo) is solution of ( [Tol ). Conversely, if (uq, w) is solution 
of (jToJ then u defined by 

u\bq = Wo 

m|b± = (/?, w, (^), where ((9 = uo|p± 



associated to the same value uj is solution of ( jg^j i. Moreover, the multiplicity of w is the same 
for the two problems. 



Whereas the problem \£p\ was linear with respect to the eigenvalue uP but defined on an 
unbounded domain, the problem dJoJ is set on a bounded domain but non linear. 



Remark 4.5. Note that the problem \Eq) is also non linear with respect to f3 (whereas 



the problem (£p 1 can be rewritten as a quadratic eigenvalue problem). In other words, this 
difficulty would be present if we decided to fix lo and look for /3. 

We now introduce the solution algorithm of the non linear eigenvalue problem, explain how 
to compute the DtN operators in the case where ^ aessiP) and finally deduce another 
characterization of the essential spectrum <Tess{l3)- 
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4.2. Solution algorithm. For c? ^ aess(/3), we denote by v4o(/?, a) the operator 

Ao(^,a) = - — A 
Po 

+a:,uo + A+(^,a)Mo = 0, onr+ 



i^(Ao(/3,a)) = {7.oeif'(A,So), 



-9:r-ito+A a) 1*0 = 0, on , }. 



Proposition 4.6. 77ze operator Ao(/?, a) is self-adjoint, bounded from below andwith 
compact resolvant. 

Proof. To show that Aq{(3^ a) is self-adjoint, the idea is to use the definition ( |4.1| i of 
the DtN operators. The only difficulty now is to prove that ^o(/3, oi) is bounded from below. 
Using again the relation ( |4.1| l, we obtain easily 

Vu e Z?(^o(^,a)), (Ao(/3,aKw) = ||Vu|li.(B„) + \\yu+\\l,^s+) + ll Vu" H^.^s-) 

\Jn+ Ja- 
where = W^iP, a; u|p+) and u~ = u~{f3, a; u|p- ) satisfy 

C C" 

- dist(a,a±(^))""'rj"«^'^(r^) ^ dist(a,a±(/3))"""^^(^») 

using the well-posedness of the problems (T'*) and the continuity of the trace application 
from H^{Bq) onto iJ^/^(rJ). It follows that the operator is bounded from below 

□ The spectrum of a) is then a pure point one and consists of a sequence of eigenvalues 
(/!„(/?, a))n of finite multiplicity tending to +oo: 

Pi{fi, a) < A^2(/3, a) < ... < Pn{fi, a) -> +oo 

The explicit expression of these eigenvalues using the Min-Max principle yields some regu- 
larity properties of each eigenvalue with respect to a. 

Theorem4.7 (Characterization of the eigenvalues (/x„(/?, a))„). Letusfixf3 e [0,TT/Ly[ 
and a ^ aessiP) cind denote Hq = L^{Bo, podxdy). In the following Vm is the set of sub- 
spaces of dimension m of 

HliBo) = {uo e H\Bo), Molso - e'^uolgj 

Then Vm G N, 

/b,J^^oP + (A+(/3,a)Mo,wo)r+ + {^^ {P^oi)uo,uo)y~ 



^„j(/3,a)= inf sup 

and < •, • >p± denotes the duality product between Hp^^^{T^) and H^^^{T^). 
Moreover, the functions a) i— >■ fJ-miP, ct) are continuous. 
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Proof. The first part of the theorem is a classical application of the min-max principle 
since Ao{P, a) is selfadjoint and bounded from below. 

For the continuity, let us first note that for any /3, it exists a constant C{a) with depends 
continuously of a such that for any a 

yu e DiAo{/3, a)), (Ao(/3, a)u, u) + C(a)||u||2,,(^^) > \\u\\l^^Ba) (4-3) 

(see the proof of Proposition |4.6| l. Let us then consider the following Rayleigh quotients : 

/bo'^^oI^ + (A+(/3,q;)uo,mo)p+ + (A-(^,a)uo,uo>r- + C {a)\\uo\\ (Bo) 



7^(/3,a;uo) = 



I "oil Ho 



For any (3, let ai and a2 be such that the DfN operators a^; •) are well defined. For aU 

uq e H^{Bq) such that ||uo||ifo = 1> we have 

|7^(/3, ai; uo) - 7^(^, uo)| < |((A+(/3, ai) - K+{P, aa)) uo, uo)r+ 1 

+ |((A-(/3, ai) - A-(/3, aa)) uo, wo)r- I + |C(ai) " C{a2)\\\u\\Hi(Bo) 

So using Proposition |4.3| and its proof and the continuity of the trace application from (Bq) 
onto iJi/2(r±), 

we obtain the existence of a constant C independent from ai and a2 such 

that 

C 

dist(Q!f , CT+(/3)) dist(a^, cr+(/3)) 

+ \C{ai) - C{a2)\\\uo\\HHBo) 

Taking into account finally the coercivity property ( |4.3| l, we obtain 

C" 

n{(3,a2;uo) < TZ{l3,ai;uo) + ^ tawA- 2 7^l"i ~ al\TZ{l3 , ai; uq) 

dist(a^, a+[p)) aist[a2,(J+(p)) 

+ |C(al)-C(a2)|7^(/3,al;^io). 

For all m, if we denote 

Vto e N; a)= inf sup 72.(/?,a;uo) 



we obtain from the last inequality for all m 

dist(a2,a+(/3))dist(a2,a+(^)) 
+ \C{ai)-C{a2)W,nil3,ai) 



Vm{l3,a2) < V,n{l3,ai) + ^.^^f^2 IR\\A- t( 2 , raw l"l ~ "2l^m(/^'"l) 



By exchanging the roles of ai and a2, we finally obtain 

^dist(a2,a+(/3))dist(a2,a+(/3)) 
+ \C{ai) - C{a2)\ max(i/,„(/3,ai),i/,„(/3,a2)) 



ttl) I < 3JT777T2 ( iQW A' t( 2 Z t aw 1*-^! ~ '-^2 1 (fm(/3, Qfi), Vmil^, 0^2)) 



We deduce the continuity of the /im(/3, a), C(q;) depending continuously of a. □ 
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Consequently, using Theorem 4.4 the solutions of the problem ( i are obtained as explained 
in the following proposition. 

Proposition 4.8. The solutions of the non-linear problem are the roots of the 
equations : 

^'^ t cFessiP) and 3to > 1, = w^. 

We then infer the iterative algorithm for the computation of the guided modes and associated 
eigenvalues with two nested loops: 

• the outer loop consists in a fixed point algorithm to solve the non hnear equations: 



/i™(/3,a;) = cj^, cj^ ^ o-ess(/3); 



• each iteration of this fixed point algorithm requires the computation of the m— th 
eigenvalue fimiP, a) of the operator Ao(/3, a) (and possibly the derivative of fim{l3, a) 
with respect to a if a Newton method is used to solve the fixed point problem). 
This algorithm is quite classical for the computation of guided modes in open waveguides 
(see 111). Here the novelty comes from the facts that the eigenvalues could belong to any 
gap of the spectrum and that the operators {/3, u) have no analytical expression. We show 
in the following section how they can be computed in practice. 

4.3. Characterization of the DtN operators . As explained in [22], the construction of 
the operators A* a) is based only on the resolution of a family of cell problems and the 
resolution of a stationary Ricatti equation. For the sake of clarity, we shall recall the method 
for the construction of A+. 

We will suppose again in this section that ^ a^ssiP) and the problem (V^) is well posed. 
We shall use the division of the half-band into periodicity cells separated by vertical seg- 



ments (See Figure 4. 1 



+ C30 



B+=[jCn, C„:=Ci + ((n-l)L„0), 



(4.4) 



The segments r„ = Fq + {nL, 0), with Fq = F+ can all be identified to the leftmost one 
Fq ~ [—Ly/2, Ly/2] and the cells C„ can all be identified to the first cell Ci = C. By 





Fig. 4.1. Notation in a half band 

periodicity in x, the construction of the unique solution u^{/3, a, (p) of (V^) will be reduced 
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to the knowledge of a linear operator, called the propagation operator, denoted -P(/3, a) and 
defined by 



P(/3,a): H'/'iTo)^ h'/\To) 



(4.5) 



One can show that P(/3, a) for any /3 and a is compact (using interior elliptic regularity for 
a, ip)), injective (using an argument of unique continuation) and has a spectral radius 
less than 1 when ^ (Tess(/3) (because of the nature of u^{f3, a, ip)). See 1221 [13 for 
more details on the proof of these results. 

Using the periodicity of the problem, one easily see that 

VneN, V(^ei?y'(ro), u+(/3,a,(p)^ =(P(/3,a))" ^. 



r„ 



1 /2 

Then by linearity, we haveV?! G N, G H^' {Tq) 

w+(/3,a,^)|^ =e°(/?,a,(P(/?,a))"-^ (p) +ei(/3,a,(F(/3,a))" ^) (4.6) 

where for all Lp in 77^^^ (Fq) eo(/3, a, tp) and ei(/3, a, (yc) are the solutions in H^{C) to the 
following elementary cell problems 

- Aee~ Ppa^ee^O in C, ^€{0,1} 

satisfying /3-quasi periodic boundary conditions on Si — Y.OC and Si = S n C 

and Dirichlet boundary conditions on and F^ 



(4.7) 



(4.8) 



eolro = (p 
eolri = 



and 



ei|po = 
eilri = V 



(4.9) 




Fig. 4.2. The periodicity cellC 

Remark 4.9. The cell problems are well posed except for a countable set of frequencies 
which correspond to the eigenvalues of the operator pp defined in a cell with Dirichlet 
boundary conditions. These forbidden frequencies are introduced because of the Dirichlet 
boundary conditions and using Robin-type cell problems would solve this issue. We choose 
here Dirichlet-type cell problems only to simplify the presentation of the method. 
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Formula ( |4.6| ) shows that the computation of the solution is achieved through the char- 
acterization of the operator P{(3,a). At this stage of the exposition, the definition of this 
operator relies on u^{P, a, ip) which is a solution of a problem posed in an unbounded do- 
main. We shall see in the following how to determine this operator by solely solving local 
problems of the type (|4.7||4.8[|4.9|l, which is one key point of the method. 



Note that the relation ( |4.6| l ensures that m+ (/3, a, is the solution of the Helmholtz equation 
inside each cell C„. To make the characterization complete, we have to add the transmission 
condition across Tj (the continuity of the normal derivative of u+ accross Tj), which gives 

Tio(/3, a)P{(i, af + [T^oifi, a) + Tn(/3, a)) P{f3, a) + Toi(/3, a) - 0. 

where the operators Tpq{/3, a) are four corresponding local DtN operators : 



Toi{(3)(p = -dxeo{/3,(p) 



po 



ri 



rio(/3)(^ = -a,ei(^,(^) 



(4.10) 



Actually, for any f3, if ^ aessiP), this equation characterizes uniquely the operator 

P(/3,a): 

Theorem 4.10 (Characteristic equation). Suppose that o? ^ (Tess(/?) and the prob- 
lem (V'^) is well posed. The operator P{l3,a) is then the unique compact operator of 
Ci^Ha^"^ (Tq)^ satisfying the condition 



p{P{P,a)) < 1 
which solves the stationary Riccati equation: 

where 



(4.11) 



(f^(/3,a)) 



X ^ Tio(/3,a)X2 + (roo(/3,a) +Tu(/3,a))X + Toi(/3,a). 

Proof. As explained before, the operator P{/3, a) is a compact operator whose spectral 
radius is strictly less than one. Moreover, the continuity of the normal de rivative of the unique 
solution ( 



, a; if) accross each Tj, P{P, a) is solution of i£ {/3, a) 



Reciprocal ly, let P be a compact operator with spectral radius p{P) < 1 which satisfies 
i£^{f3, a) I. The function defined by 

u{^)\^ = e" a, ^) + el {f3, a, (Pf ^) 

satisfies the helmoltz equation cell by cell, is continuous by definition a nd has con tinuous nor- 
mal derivatives accross each Tj because the operator P is solution of {£^{(3, a) i. Moreover, 
the fact that p{R) < 1 implies that u{(p) belongs to H^{B~^). Indeed, the property 



lim II P 



.ri II 1/n 

£(Hi/2(r)) 



P{P) 
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implies that, for some p* G {p{P), 1) and N large enough 

Vn > N, WniciHyHD) < P* 
so that, it exists C independant of N, such that 

We have then constructed a solution in (A, 5+) of (V^). By well-posedness of this prob- 
lem, the solution is unique so u = u+ and then P = P{f3, a). □ 

Once P{/3,a) is determined solving the stationary Ricatti equation, we build cell by cell 
the solution u+ using ( |4.6[ ) and finally using again ( |4.6[ ) for n — 0, we see that 

A+(^, a) =. Too(^, a) + T^W, a) P{f3, a) (4.12) 

4.4. Characterization of the essential spectrum. The characterization of the essential 
spectrum is a by-product of the determination of the DtN operators. 

Indeed, in [22J, we have shown that when lies in the essential spectrum (Tess(/?), not 
only the spectral radius of the propagation operator is equal to 1 but in addition this operator 
is no more the unique solution of the Ricatti equation (fij). Consequently, we can give this 
characterization of the essential spectrum of A. 



Theorem 4.11. The Ricatti equation (£r) has a unique solution whose spectral radius 
is strictly less than one if and only if o? (ji aessil^)- 



Proof. In Theorem 4.10 we have already proven that if a ^ fJessW), the Ricatti equa- 
tion (Eli) has a unique solution whose spectral radius is strictly less than one. 

For e CTess(/3), we can show that it does not exist solution of (En) with spectral ra- 
dius strictly less than one but solution whose spectral radius equal to one. Moreover, there is 
no more uniqueness of this kind of solutions. Indeed, by definition of the essential spectrum 
(see the proof of Proposition [3T|l 

i <Jessi0) ^ 3ko e] - n/L^, tt/L^] = fco) 

where (/3, fc) is the q—th eigenvalue of the operator Ap{f3,k) defined in the proof of Propo- 
sition 3.1 Let us denote eq{/3, k) a corresponding eigenvector. Using the same idea than in 
the proof of Proposition 4.3 we could show that k ^ ojq{/3, k) is an even function and then 

ko) = UJq{f3, -fco) = 

Let us denote now (p* — eq{f3, ±fco). Using the well-posedness of the cell problems, the 
eigenvectors eq{/3, ±fco) can be expressed by 

e,(/?, ±fco) = e"(/3, a, ^±) + e±^'=°^== e\/3, a, ^±) 

and, thanks to their quasi-periodicity, can be extended as H^(A, 5+) function by 



eg(/3, ±fco) 



.M"'-^'^'"'^'^ e°(/3, a, ip^) + e±^"'=«^- e\p, a, <p^). 
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We have found then solutions of the problem {V^) with (p* as Dirichlet Boundary condition 
on Tq. We can then construct at least three solutions of the Ricatti equation whose spectral 
radius is equal to one : 

Pi such that Pi ip+ = Q^koL^ ^+ 

P2 such that P2 ^p' 



P3 such that P3 ip+ = e*''o-^=" (p+ and P3 ip^ = 



□ 



5. Numerical results. We want to apply this approach to compute the guided modes of 



the media, whose refraction index is given by (see the isovalues Figure 5.1 

y{x,y) e [-0.5,0.5]2,Vn, to€ Z*, 

The period here is equal to 1. 



Pp{x + n,y + m) = 1 + 16exp( ^ ) 



Po{x,y + m) = 1 



0.22 



(5.1) 



>.4 




Fig. 5.1. The lineic perturbed periodic media. 



For any j3 e [Q^n/Ly] and e [0,20], we determine if is in CTess(/3) or not. If not, 
we compute the Dirichlet-to-Neumann operators A^(/3, a). The corresponding algorithm is 
the following : for each /3 and a 

1. solve the cell problems (| 4.7|4.8|4.9 1 - using a FE method or a mixed finite ele- 
ment one (see [2?, "TS] [Tif for more details) and compute the local DtN operators 
T±(/3,a) (see CT ); 

2. solve the stationary equation ( |4.10| l using a spectral decomposition method (which 
leads to a quadratic eigenvalue problem) or a modified Newton algorithm (see 1221 
[15] El for more details): 

(i) If it exists a unique solution whose spectral radius is strictly less than one then 

^ <Jess{P), SO we can determine the propagation operator P*(/3,a) and 
deduce the DtN operators A='= (/?, a) by ( |4.12| i. 

(ii) If it exists a solution whose spectral radius is equal to one then G <Jess{P) 
and we stop the computation. 
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We can then, for any f3 and a, compute the eigenvalues Hm{P, ct). 



We plot on Figure 5.2 the isovalue lines of the function log|/ii(/3, a) — a \ for (3 £ [0, tt/Xj, 
anda2 g [Q, 20] \ ae«s(/3). 




0.5 1 1.5 2 2,5 3 



y3 

Fig. 5.2. Isovalues o/log|/xi ((.j) — aj^|; /9 £ [0, tt/L], G [0, 20]. The two Mack circles correspond to 
guided modes represented F/gMre [Tj] 

For a fixed /?, the white regions corresponds to the essential spectrum aessiP)- A part of the 
dispersion curves are given by the blue lines corresponding to the null isovalue /ii(/3,a;) — 
ijp — 0. The main difference with homogeneous open waveguides is that the function 
fimiP, w) — may vanish several times for a fixed /3. 



We want know to compute for fixed /3, the eigenvalues and the corresponding guided modes. 
To do so, we use the following steps : 

1 . Solution of the non linear eigenvalue problem : using the algorithm of the non 

is computed so that 



4.2 



linear eigenvalue problem described at the end of Section 
it exists m such that uj) — uP' and uj^ ^ aessiP)- Let uq be an eigenvector of 

Ao{/3, ui) associated to the eigenvalue fimiP, 
2. Construction of the guided modes in one band : by continuity arguments, the 
restriction to a band i? of a guided mode u is given by 



\B + 



uo 



= oj, with (^9+ 

= u~{/3,uj,(p~) with Lp~ 



uq\ 



uo\ 



where for any given (p G H^^^(T^), u^{/3,uj,ip) is the unique solution of ( V^ i. To 
compute these solutions, we have to 
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• solve the cell problems ( 4.7|4.8|4~9 1 for a — u and compute the local DtN 
operators T^q{P,uj) (see ( |4. 10| i); 

• solve the stationary equation ( |4. 10| i. Since lo^ ^ a^ssiP), it exists a unique 
solution whose spectral radius is strictly less than one which is the propagation 
operator {(3, lo) 

• reconstruct the solutions cell by cell using \A.6\ . 

3. Construction of the guided modes in the whole domain : by /3— quasi-periodicity, 
the guided modes can be constructed "band by band" : 

V(x,y) e B, VqeZ, u{x,y + qLy) = u{x,y)e"'^^y. 



We then represent Figure |5.3| guided modes for two values of /3 in eight periods from each 
side of the line defect. Figure 5.3(a) corresponds to the eigenvalue in the first gap of A{0.5) 



with a well confined mode whereas in Figure 5.3(b) the eigenvalue belongs to the fourth gap 
of A{1.9)) and the associated guided mode is not well confined. 



Ill mil ^ 
lull ^ 



(a) 13 = 0.5, a;2 = 3.465 



'ifSfifll M M II 1 II M «i M «i II II 
l< I* I* l> •* Il II I II II II «l il •! •( 
• I Il II Mil II I nil II II II ti I* 
• • !• Iili^il •! il •! I 



2 4 



(b) 13 = 1.42, ^2 = 10.46 
Fig. 5.3. Modes: well confined (left); not well confined (right). 



6. Conclusions. We have shown that this method is well adapted to the computation of 
guided modes which are not well confined. The advantage of this method is that it is exact 
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independently from the confinement of the computed guided modes. It can be well adapted 
for the study of the dispersive curves near the essential spectrum. In terms of computation 
time, the solution of the non linear eigenvalue problem is much more costly than the one of a 
linear eigenvalue problem. However, on the contrary of existing methods, we do not have to 
compute independently -and often using a separate code- the essential spectrum of the media 
: here it is a by-product of this method. 

In a forthcoming article, we will compare exactly this method with the Supercell Method 
in terms of accuracy and computation time. Moreover, this DtN strategy will be apphed to a 
;9-formulation which is more adapted for dispersive media. 

Another interesting perspective concerns the study of the dispersive curves : with our charac- 
terization, we could study if they can or cannot be constant with respect to P (and conclude 
on the existence or not of boimd state) and analyze their behavior near the essential spectrum. 
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